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From Weak Lensing to non-Gaussianity via Minkowski Functionals 
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ABSTRACT 

We present a new harmonic -domain approach for extracting morphological information, in the 
form of Minkowski Functionals (MFs), from weak lensing (WL) convergence maps. Using a per- 
turbative expansion of the MFs, which is expected to be valid for the range of angular scales probed 
by most current weak-lensing surveys, we show that the study of three generalized skewness pa- 
rameters is equivalent to the study of the three MFs defined in two dimensions. We then extend 
these skewness parameters to three associated skew-spectra which carry more information about 
the convergence bispectmm than their one-point counterparts. We discuss various issues such as 
noise and incomplete sky coverage in the context of estimation of these skew-spectra from realis- 
tic data. Our technique provides an alternative to the pixel-space approaches typically used in the 
estimation of MFs, and it can be particularly useful in the presence of masks with non-trivial topol- 
ogy. Analytical modeling of weak lensing statistics relies on an accurate modeling of the statistics 
of underlying density distribution. We apply three different formalisms to model the underlying 
dark-matter bispectrum: the hierarchical ansatz, halo model and a fitting function based on numer- 
ical simulations; MFs resulting from each of these formalisms are computed and compared. We 
investigate the extent to witch late-time gravity-induced non-Gaussianity (to which weak lensing 
is primarily sensitive) can be separated from primordial non-Gaussianity and how this separation 
depends on source redshift and angular scale. 
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1 INTRODUCTION 



Since the first me asurements were published dBeacon. Refregier&EUisI l200d : IWittman et ail |2000| : _K aiser. W ilson & Luppinol 120001 : 
I Waerbeke et all l2000) there has been tremendous progress in the field of weak gravitional lensing, regarding analytical modelling, as well as 
technical specification and control of systematics in observational surveys. Ongoing and planned weak lensing surveys fsee lMunshi et al.l d2008h 
for a review) such as the CFHlF] Legacy Survey, Pan-STARR50, the Dark Energy SurvejQ and, further in the future, the Large Synoptic Survey 
Telescopfl JDE1V0 and EuclidQ will map the cosmological distribution of dark matter and probe the properties of dark energy in unprecedented 
detail. Owing to the greater sky coverage, tighter control on systematics and increased number-density of source galaxies it will b e soon possible 
to extract higher-order statistics (i.e. beyond the two-point correlation function), such as multispectra; see e.g. |Pen et al. I J2003I) . Non-linearity 
induced by gravitational effects is generally used to break the degeneracy between the amplitude of matter power spectrum as and the matter 



1 http://www.cfht.hawai.edu/Sciences/CFHLS/ 

2 http://pan-starrs.ifa.hawai.edu/ 

3 https://www.darkenergysurvey.org/ 

4 http://www.lsst.org/llst_home.shtml 

5 http://jdem.gsfc.nasa.gov/ 

6 http:/, 'http://sci.esa.int/euclid 
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Figure 1. The power spectrum is plotted as a function of the harmonics I in the left panel. A WMAP7 background cosmology is used. The results are displayed for 
three different source redshifts z s = 1.5, 1.0 and 0.5. The cosmological parameters are T = 0.1956, Q = 0.279, A = 0.721 and erg = 0.817 respectively. The 
dot-dashed curve shows predictions from the halo model for the same WMAP7 background cosmology and for the source redshift z s = I. The diagonal entries of 
the bi-spectrum are plotted as a function of the harmonics I in the right panel. The results are for z s = 1.0. Two different approaches are persued in computation 
of the bispectrum. The bispectrum results based on extnsions of perturbation theory Eg. ( 1341 are plotted using solid lines, the halo model predictions i.e. Eq. j3U are 
shown using dashed lines. 



densi ty parameter Qm; three-point statisti cs such as the bispectrum (the three-point multispectrum) are the best studied statistics for this pur- 
pose (Villu msenlll 9961 .1 Jain & Seliakll 19971) . Weak lensing can therefore play an important role in breaking degeneracies, which makes it an ideal 
complement to Cosmic Microwave Background (CMB) studies and studies involving large scale structure (LSS) surveys. 

Two-point statistics, principally the power spectrum, of density perturbations remain the most frequently used statistical tool for many 
cosmological studies. Weak lensing surveys probe the non-linear regime and are therefore sensitive to non-Gaussian signatures which can not 
be probed using two-point statistics . The statistics of shear or convergence probe the statistics of underlying ma s s distribution in an unbi- 
ased way ijj ain, Selia k & White]|200d:lMunshi & Jaidl200ll;lMunsIrfe00rj|;|Munshi & Jaidl2(Ml200ll:IValapadl200d : lMunshi & Valageasll2005l : 
IValageas. Barber. & Munsh i 2004; Valage as. Munshi & Ba rber 2005; Taka da & Whitell2003l ; iTakada & Jainll2004t) . but are very sensitive to non- 
linear evolution driven by gravitational clusterin g. A number of analyti cal sche mes, from perturbative ca lculations to halo models have therefore 
been employed to model weak lensing statisti c s iFrvl (119841); ISchaefferl j 198 4); Bernardeau & Schaeffer ( 1992); Szapudi & Sza layl d 19931 II 997); 
IMunshi et aj dl999h : iMunshi. Coles & MelotJ ll999ij |bl); iMunshi. Melott & Colesj dl999h ; IMunshi & Coles! ( feOOd . I2002L 120031) ; Icoorav & Sethi 
d2002l) ). In addition to studying the statistics in projection on the sky, they have also been studied in three dimensions using photometric red- 
shifts. It has been demonstrated t hat this approach can tighten observational constraints on such quantities as the neutrino mass and the dark 
energy equation of s tate parameter dHeavensll2003l ; lHeavens et alj|200p| ; iHeavens et alll2006l ; lHeavens. Kitching & Verde 2 0071 ; ICastro et alll2005l ; 
Kitching et al. 2008). To mographic techniques have also be en employed as an intermediate strategy between projected surveys and 3D mapping 
dHuj|l999l ; lTakada & Jainll2004l2003l ; ISemboloni etal2008h . 

Minkowski Functionals (MFs) are morphological descriptors that are commonly used in many cosmological cont exts. They can be 



Komatsu et al. 2003; Eriksen et al. 2004), weak lensing surveys (Matsubara & Jain 2001; Sato et al. 2001; Taruya et al. 2002) and galaxy surveys 


jGott. Melott & Dickinson! 1 19861: IColesjl 19881: iGott et alj 19891; iMelott et alJI 19891: iMoore et alJ 


1992; lGottetal.ll 19921: iRhodas. Gott & PostmarJ 


Il994l:lcanavezes et al. 1998; Park et al. 2005; Hikaae et al. 2008. 2006; Hikage. Taruva & Sutol 


2003l:lffikage et alj|2002b. Unlike the multispec- 



tra, discussed above the topological descriptors carry information of all orders (in a statistical sense). In the context of CMB studies, the MFs are 
used to probe primordial non-Gaussianity. For large scale structure studies using projected or redshift galaxy surveys, the non-Gaussianity probed 
is mainly that which is induced by gravity. While galaxy surveys suffer from uncertainties relating to the nature of galaxy bias, weak lensing 
surveys will provide an unbiased probe to probe the clustering of dark matter. The MFs will be an important tool in this direction, along with other 
statistics that can be used to probe non-Gaussianity to break the parameter-space that are unavoidable when the power spectrum alone is used. 

This paper is organized as follows. In §2 we review the formalism of Minkowski Functionals. In §3 we link the statistics of weak lensing 
convergence and the underlying density distribution. In §4 we introduce the concept of generalized skew-spectra and show how these power 
spectra can be used to study the Minkowski Functionals. In §5 we review the analytical models that are typically used for modelling of dark matter 
clustering. 
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Figure 2. The hermite polynomials that are used as a basis function for expanding the Minkowski functional in the weakly non-gaussian limit. The plots show 
Hi (x), (solid) H2 (x), (short-dashed) ^3 (x), (long-dashed) H4(x)(dot-dashed) as a function of the argument x as depicted. 



2 FORMALISM 

The MFs are well-known morphological descriptors which are used in the study of random fluctuat i on fie lds fields. Morphological properties are 
defined to be those properties that remain invariant under rotation and translation; see lHadwigerl l fl959l) for a more formal introduction. They 
are defined over an excursi on set S for a give n threshold v. The three MFs that are defined for two dimensional (2D) studies can be expressed, 
following the notations of (Hikage et al. 200_8|), as: 



Va(v) = f da; Vi(v) = \ [ dl; V 2 (v) = ^- [ Kdl (1) 
Je 4 J as 27V J 812 

Here da, dl are surface are and line elements for the excursion set E and its boundary <9E respectively. The MFs V k (y) correspond to the area of 
the excursion set E, the length of its boundary <9E as well as the integral of curvature K, along its boundary which is also related to the genus g 
and hence the Euler characteristics \- 

In our analysis we will consider a smoothed random field k(Ci), with mean {k(Q,)} — and variance <r 2 = ( k2 (^)}; f° r me ti me being 
k is a generic 2D weakly non-Gaussian random field defined on the sky although we will introduce more specific examples later on. The spher- 
ical harmonic decomposition, using Yj m (fl) as basis functions, k(Ci) — Xw m K imYi m (Cl), can be used to define the power spectrum Ci using 
(KimKy^i) = Ci&ui8 mm i which is a sufficient statistical characterization of a Gaussian field. For a non-Gaussian field, higher-order statistics such 
as the bi- or tri-spectrum can describe the resulting mode-mode coupling. An alternative to this laborious expansion in multispectra, topological 
measures such as the Minkowski functionals can be employed to quantify deviations from Gaussianity and it can be shown that the information 
content in both descriptions is the same. At leading order the MFs can be constructed completely from the knowledge of the bispectrum alone. We 
will be studying the MFs defined over the surface of the celestial sphere but equivalent results can be obtained in 3D using Fourier decomp osition 
(Mun shi 2010, in preparation). The behaviour of the MFs for a random Gaussian field is well known and is given by Tomita's formula jTomital 
1986). The MFs are denoted by Vk{y) for a threshold v — k/uo, where ctq = (k 2 ) can be decomposed into two different contributions, Gaussian 
(V fe G (^)) and non-Gaussian (<5I4(V)), i.e. Vk{v) = V k G (u) + 8Vk(v). From our perspective we will be more interested in the non-Gaussian 
contribution, i.e. 6Vk (^)- We will further separate out an ampli tude A in the expres sions of both of these contributions which depend only on the 
power spectrum of the perturbation through ctq and o\ (see e.g. lHikage et alj|2008h ): 



V k G {v) = Ae* V 



Hk-i; 



&V k {v) = Aexp (~\ [<SVf V^o + 5Vt 3) {v)al + 8V k ii] \u)ol + •••] 



ls m H k+2 (u) + ^S w H k (u) + 



6 



A = 



( 27r )(fc+D/2 u2-ku k \V2a 



(2) 
(3) 



The constant ui k introduced above is the volume of the unit sphere in k dimensions, i.e. w k = n k ^ 2 /T(k/2 + 1); in 2D we will only need 
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wq = 1, u)i = 2 and ui2 = 7T. The lowest-order Hermite polynoimals Hk(y) are listed below. As mentioned previously, the expressions consist 
of two distinct contributions. The part which does not depend on the three different skewness paramters S^°\ signifies the MFs for 

a Gaussian random field. The other contribution 5Vk(y) represents the departure from the Gaussian statistics and depends on the generalised 
skewness parameters defined in Eq.((5} and Eq.((9). We have expanded the total non-Gaussian contribution into perturbation series in uq. While the 
lowest order terms <5 V fc (v) are determined by various one-point moments related to the bispectrum the next-to- leading-order terms <5 V fc (v) are 
connected to similar one-point moments related to trispectrum (also known as the kurtosis). In projected surveys, even for relatively small angular 
smoothing scales, the leading order terms are sufficient to describe the non-Gaussian departures in the smoothed convergence field k(8 3 ). 



ff " lM = V? e ' 2/2erfc (^) ; HoW = 1 ' H ^ v ) = v > 

H 2 {v) = v 2 -l, H 3 (u) = u 3 -3u, H 4 (u) = v A - 6u 2 +3. 

^M = (-irexp(^)-Ae X p(^) (4) 

The various quantities Oj that appear in Eq.((3} can be expressed in terms of the power spectra Ci and the shape of the observational beam 6j, 
The moment Co is a special case for which Oq corresponds to the variance. The quantities oi, 02 are natural generalisations of variance, with 
increasing j corresponding to increase weight towards higher harmonics tr| = ^2(21 + + l)] J Ci& 2 . The variance that will mostly be used 
are a 2 = (k 2 ) and <y\ = <(Vk) 2 ). 

The real-space expressions for the triplets of skewness are given below. These are natural generalisations of the ordinary skewness 5° 
that is used in many cosmological studies. They all are cubic statistics but are constructed from different cubic combinations: 



The expressions in the harmonic domain are more useful in the context of CMB studies where we will be recovering them from a masked sky 
using analytical tools that are commonly used for power spectrum analysis. The skewness parameter S^ 1 ' is constructed from the product field 
k 2 and V 2 /, whereas skewness parameter 5 1 ' 2 ' relies on the construction of [Vk • V/t] and V 2 k. By construction, the skewness parameter 
has the highest weight at high I modes and has th e lowest weights on high I modes. The expressions in terms of the bispectrum i?i 1 i 2 i 3 (see 
Eq.lllOt for defintion) take the following form (see e.g. lHikage et alj J2OO8I) ): 



S< ~ K ' = ■^^2 B hhh I lil2h (6) 

k 

S (»»v»») = ^ [ h ( h + 1) + j a (i a + 1) + l 3 (l 3 + l)]B lll2h I hhh (7) 

k 

s (v K -v.v 2 K ) =±-^2 [piGi + 1) + *a(i2 + 1) - Z 3 (k + 1)]I 3 (I 3 + 1) + cyc.porm.]iB ll j 2 j3/, 1 , !l ,3W 1 Wi a Wi a (8) 



^ (2Zi + 1)(2Z 2 + 1)(2Z 3 + 1) / h h h 
lllhh V 4tt I 



(9) 



The bispectrum Bi t i 2 i 3 used here defines the three-point correlation function in the harmonic domain. A reduced bispectrum &i 1 i 2 ; 3 can also be 
defined which can directly be linked to the flat-sky expressions. 



(Ki imi Ki 2m2 Ki 3 m 3 )c = [ mi m2 mjj J Bhhh'i Bhhh = hiiihbiiiih- ( 10 ) 

The expressions for the MFs in Eq.([3j depend on the one-p oint cumulants . Howe ver it is possible to define power spectra associated with 
each of these skewnesses following a procedure developed in iMunshi & Heav ens (2010). This will mean we can also associate a power spectrum 
with V^ 3 ' which will generalize the concept of MFs in a scale-dependent way. The power spectrum that we associate with MFs will have the same 
correspondence with various skew-spectra S; as the MFs have with one-point cumulants or S 1 ' '. The power spectra so defined will however 
have more power to distinguish various models of non-Gaussianity. This is one of the main motivations behind generalising the concept of MFs, 
each of which is a number, to a power spectrum, which contains scale information. 

The series expansion for the MFs can be extended beyond the level of the bispectrum; the next-to-leading-order corrections terms are related 
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to trispectra of the original fields and various derivatives constructed from them using differential operations such as V • V, V 2 . These corrections 
are expected to be sub-dominant in the context of CMB studies for the entire range of angular scales being probed. 

The results here correspond to analysis of convergence maps, which are spin-0 objects. It is possible to extend these results to spin-2 fields. 
Such results will be interesting for analysis of weak lensing shear and flexions, but a detailed analysis will be presented elsewhere. 



3 CONVERGENCE POWER SPECTRUM AND BISPECTRUM 

The convergence n(Cl, r s ) can be treated as a line-of-sight projection of the density contrast S(£l, r) along the direction Cl (r is comoving radial 
distance) out to a source redshift z s (r s ) with a redshift-dependent weight function uj(r, r s ): 



r r * - 3 Hn d,A(r — r s ) 

K(Q,r s ) = drw(r,r a )5(Q,r); U)(r,r s ) = - , s ■ r s = min(ri, r 2 ). (11) 

Jo 2a c dA(r)d,A{r s ) 

The weight functions cj(r) for weak lensing depend on the angular diameter distance dA(r), Hubble constant Ho, matter density parameter Qm 
and the scale factor of the Universe a = 1/(1 + z) at a redshift z. The angular diameter distance dA(r) is linked to the total matter content 
fi and the Hubble constant H , i.e. d A (r) = KT 1/2 sin(KT 1/2 r), KT 1/2 sinh((-K)- 1/2 r), r for open, closed and flat Universes; here 
K = (fio — 1)Hq. We will consider the projected cross-power spectra Ci that depend on two differ ent redshift z\ and 22 which is a f unction of 
the underlying matter power spectra Ps(k, r) which, in the small-angle approximation (Limber 1954), can be expressed as jKaisej|l992h : 



Cl= L dr ^r p *{dA^y r )-> 

Analytical modelling of the convergence bispectrum B depends on modelling of the underlying matter bispectrum Bs : 



Jo d A( r ) \ d A(r) d A (r) d A {r)J 

we will discuss the analytical models we use to construct Bg in later sections. This equation can also be used to express the reduced bispectrum 
bhhh introduced before in Eq. J 10t . Estimation of individual modes of the bispectrum defined by specific choice of the triplets (h, h, h) is difficult 
when the data are noisy, but it is possible to extract the cross-correlation of product maps k 2 (Q.) against k(CI). If we denote the harmonics of the 
product map K 2 (Cl) as /-c^ = J dCl Y* m (Cl)K 2 (Cl) and similarly K lm = J dCl Y^(f2)/t(f2); then the associated power spectrum is constructed as 

^? ~ 2IT1 ^ e ( K Wn K lm) ' s called the skew-spectrum dCooravlboOll) . We will next generalize the concept of skew-spectrum and introduce 
a set of generalized skew-spectrum that can be used to construct the Minkowski Functionals at the lowest level of non-Gaussianity. 



4 THE TRIPLETS OF SKEW-SPECTRA AND LOWEST-ORDER CORRECTIONS TO GAUSSIAN MFS 

The skew-spectra are cubic statistics that are constructed by cross-correlating two different fields. One of the field used is a composite field typically 
a product of two maps either in its original form or constructed by means of relevant differential operations. The second field will typically be 
a single field but may be constructed by applying various differential operators. All three skewnesses contribute to the three MFs that we will 

consider in 2D. 

The first of the skew-spectra was studied by dCooray|200ll) and was later generalized bv lMunshi & Heavens! j2O10h and is related to commonly 
used skewness. The skewness in this case is constructed by cross-correlating the squared map [k 2 (S1)] with the original map [«($/)]. The second 
skew-spectrum is constructed by cross-correlating the squared map [ft 2 (f2)] against [V 2 ft(f2)]. Analogously, the third skew-spectrum represents 
the cross-spectra that can be constructed using [Vre(n) • Vre(n)] and [V 2 k(^)] maps. 
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^! s ' k2 ' V2k) " I^!^TT^ Real([v2Kl!m[K2]L) 

T~2 £ ['(' + 1) + + 1) + Ia(k + 1)W; 2 JlhljWiW^Wi,, (15) 

CT (7 1 l L J 

1 1 m 

= £ [[*(* + 1) + + 1) - + l)]Ia(Ja + 1) + cyc.perm.] /3„ li2 J lllh W l W h W h (16) 



167T(To 
1 

167TCT 2 , 



h^s _ /(2*i + l)(2Z 2 + l) / h l 2 l 3 



(17) 



111213 ~ 2l 3 + 1 y (2Z 3 + 1)4tt V 

s (i) = J2( 21 + I ) s l i) 



^ = iE( 2i + 1 )[ i ( i+1 )] ic « w '' a (19) 
1 

This set of equations constitutes one of the main results of this paper. The matrices here denote the Wigner-3j symbols, Wi represents the 
smoothing window, e.g. a top hat, Gaussian or some form of compensated filter. Each of these spectra probes the same bispectrum B Hl i 2 with 
different weights for individual triplets of modes that specify the bispectrum (I, hh) and define a triangle in the harmonic domain. The skew 
spectra is summed over all possible configurations of the bispectrum keeping one of its sides at a fixed /. For each individual choice of I we can 
compute the skew-spectrum relatively straightforwardly, by constructing the relevant maps in real space (either by algebraic or differential 
operation) and then cross-correlating them in the multipole domain. Issues related to mask and noise will be dealt with in later sections, where we 
will show that, even in the presence of a mask, the computed skew spectra can be inverted to give a unbiased estimate of all-sky skew-spectra. 
Presence of noise will only affect the scatter. We have explicitly displayed the experimental beam bi in all our expressions. 

To derive the above expressions, we first express the spherical harmonic expansion of the fields [V 2 k(S1)], [Vk(O) • Vn(£l)] and [ft 2 (fi)] in 
terms of the harmonics of the original fields K; m . These expressions involve the 3j functions as well as factors that depend on various U dependent 
weight factors. 



h h I 

mi 7712 — m 



[v 2 K (f2)] ;m = J d£iY; m (Ci)[v 2 K(n)] = -i{i + i)> t ir n 

[V«(fi) • VK(£l)]im = J dnY l * m (£l)[VK(£l)-VK(£l)] = ^K hrni Ki 2 r n2 j dCl Yi* m {fl) [VY lim2 {Cl) ■ VY hrn2 (Q)} (20) 
= JX)[ii('i + l) + 'a('2 + l)-i(i + l)] / d{lYf m ({l)Y hmi ({l)Y hm2 ({l) 

liirii 

= lY / (- i riHh + l)+Hl2 + l)-l(l + l)}Khr ni Kl 2 rn 2 I h l 2 l( ^ ^ _'Y (21) 

We can define the power spectrum associated with the MFs through the following third order expression: 

Vf ' = £[V*M 2 ' + 1) = I J2(2l + 1) {s™H k (v) + k -S^H k ^{u) + fci) S ^H k . a (u) + ■■■}. (22) 

The three skewnesses thus define triplets of Minkowski Functionals. At the level of two-point statistics, in the harmonic domain, we have three 
power-spectra associated with MF V"^ 3 ' that depend on the three skew-spectra we have defined. We will show later in this paper that the fourth 
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order correction terms too have a similar form with an additional monopole contribution that can be computed from the lower order one-point 
terms in a similar way as the three skewness defined here. The result presented here is important and implies that we can study the contributions to 
each of the MFs Vk(y) as a function of harmonic mode I. This is especially significant result as various form of non-Gaussianity will have different 
/ dependence and so can potentially be distinguished from each other using this approach. The ordinary MFs add contributions from all individual 
/ modes and hence have less power in differentiating various contributing sources of non-Gaussianity. This is one of main motivations to extend 
the concept of MFs (single numbers) to one-dimensional objects similar to power spectrum. 

It is worth pointing out that the skewness and generalized skewness parameters are relatively insensitive to the background cosmology 
but quite sensitive to the underlying model of non-Gaussianity. The main dependence on cosmology typically results from the normalization 
coefficients such as ao and <j\ which are determined the power spectrum of the convergence k. 

In real space the skew-spectra can be defined through these correlation functions: 

5 (0) (ni,n 2 ) = (fc 2 (fii))Ac(n 2 )); S (1) (ni,n 2 ) = <K 2 (fii)V 2 K(fi 2 )); S (2) (ni,n 2 ) = (V«(fri) ■ Vk(Ai)V 2 k(« 2 )> (23) 
Although we have adopted a harmonic approach these correlations can equivalently be used to probe MFs especially for smaller surveys. 



5 MODELLING THE PRIMORDIAL AND GRAVITY-INDUCED BISPECTRUM 

s It is clear that we need accurate analytical modeling of dark matter clustering for prediction of weak lensing statistics, but in general there is 
no definitive analytical theory for handling gravitational clustering in the highly nonlinear regime. On larger scales, where the density field is 
only weakly nonlinear, perturbative treatments are known to be valid. For a phenome nological statistical description of dark matter clustering in 
collapsed objects on nonlinear scales, typically the halo model dCoorav & Sethll2002l) is used. We will be using the Halo Model in our study, but 
an alternative to the Halo Model approach on small scales is to employ various ansatze which trace their origin to field-theoretic techniques. Here 
we provide a quick summary of some of the analytical prescriptions that can be used to model non-linear clustering. We will also provide a brief 
description of various models of primordial non-Gaussianity arising from variants of the inflationary universe scenario. 

5.1 Hierarchical ansatz 

The hierarchical ansatz has been used for many weak lensing related work, where the higher-order correlation functions are constructed from the 
two-point correlation functions. Assuming a tree model for the matter correlation hierarchy (typically used in the highly non-linear regime) one 
can write the most general case, the N point correlation functi on, (8(r%) . . . S(r n )) c = £%(ri, . . . , r n ) as a product of two-point correlation 
functions (5(r^)<$(r 3 4 )) c = £2 (I 1 "* — r j I) l lBernardeau et all2 002). Equivalently, in the Fourier domain, the multispectra can be written as products 
of the matter power spectrum Pun(fa). The temporal dependence is implicit here. 

(JV-l) 

e JV (r lj ...,r n ) = (J(n)...c5(r n )} c = ^ Q N , a ^ ]J Uln-r^). (24) 

a, N— trees labcllings cdgcs(i.j) 

It is very interesting to note that a similar hierarchy develops in the quasi-linear regime at tree-level in the limiting case of vanishing variance, 
except that the hierarchical amplitudes become shape-dependent in such a case. These kernels are also used to relate the halo-halo correlation 
hierarchy to the underlying mass correlation hierarchy. Nevertheless there are indications from numerical simulations t hat these amplitudes be- 
come configuration-in depe ndent again as has bee n shown by high resolution studies for the lowest order case Q3 — Q ( Scoccimarr o et alll 1998b 
iBernardeau et al 2002). See Waerbeke et al (2001) for related discussion about use of perturbation theory results in intermediate scales. In Fourier 
space, however, such an ansatz means that the entire hierarchy of multi-spectra can be written in terms of sums of products of power spectra 
with different amplitudes Qn,c etc. The power spectrum is defined through (<5(ki)<5(k2)) c = (27r) 3 <53D(ki2)P*j(A;i) . Similarly, the bispectram 
and trispectrum are defined through the following expressions (<5(k 1 )f5(k 2 )5(k 3 )) c = (27r) 3 r5 3D (k 12 3)B 5 (ki, k 2 , k 3 ) and (J(ki) • ■ • <5(k 4 ))e = 
(27r) 3 <53D(ki 2 34)T' 5 (ki, k 2 , k3, k 4 ). The subscript c here represents the connected part of the spectrum and k^...^ = k^ + ■ • ■ + ki„. The 
Dirac delta functions 830 ensure translation invariance at each vertex representing the multi-spectrum. 

B*(k 1 ,k 2 ,k 3 ) Eki=0 = Q 3 [P^(fci)P* i (/% 2 ) + P„Vfci)PnVfc3) + ^ i (fc 2 )P 5 (fc3)] (25) 
T 5 (ki,k 2 ,k3,k4)2 kl= o = R a [PL{ki)Pii{k2)Pu n {k- A ) + cyc.perm.] + P 6 [P?-(fci)P?i(|ki 2 |)P*;(|ki 2 3|) + cyc.perm.]. (26) 

Different hierarchical models differ in the way numerical values are allocated to the different amplitudes. IBernardeau & SchaeffeJ dl992h consid- 
ered "snake", "hybrid" and "star" diagrams with differing amplitudes at various order. A new "star" appears at each order, higher-order "snakes" 
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Figure 3. The moments <r 2 (9 3 ) and cr\ (0 S ) (left panel) and the skewness paramters S'°) (9s), (9 S ) and (0 S ) (right panel) are plotted for source redshift 
z 3 = 1 as a function of smoothing angular scale 8 S see Eq. fll9t for definitions of cr^ and S 1 '™'. The underlying cosmology is assumed to be that of WMAP7. A top 
hat window is assumed for both of this plot. The resolution is fixed at l m ax = 4000. The underlying modelling of the convergence bispectrum Bi l j 2 ;„ depends on 
modelling of matter bispectrum Bi t j 2 j . The specific model for the underlying model that was used for this plot is based on perturbative results and its extrapolation 
to highly non-linear regime; see text for more details. The skewness paramaters are defined in Eq.((9]. The parameters Uj are defined in Eq.l ll9t . 




Figure 4. Same as previous figure but for redshift z a = 0.5 and z s = 1.5 as indicated. The skewness paramaters are defined in Eq.((9j. The parameters cr| are 
defined in Eq. )19t . For a given angular smoothing the skewness parameters increase with redshift. The cosmological parameters of the background cosmology is 
that of WMAP7. The variance parameters cr 2 increase with redshift however the skewness parameters show an increasing trend. 



or "h ybrid" diagrams are built from lower-order "star" diagrams. In models where we only have only star diagrams (Vala geas. Barber. & M unshi 
2004), the expressions for the trispectrum t akes the following form: T l5 (ki, k2, k3, ki)^ k . =0 = Q<t[P l ' (fci)P a (k2)P i (ka) + cyc.perm.]. Fol- 
lowing Valageas, Bar ber. & Munshil |2004) we will call these models "stellar models". Indeed it is also possible to use perturbative calculations 
which are however valid only at large scales. While we still do not have an exact description of the non-linear clustering of a self-gravitating 
medium in a cosmological scenario, these approaches do capture some of the salient features of gravitationa l clustering in the highly non-linea r 
regime and have been tested extensively against numerical simulation in 2D statistics of co nvergence of shear dValageas. Barber. & Munshi 2004). 
These models have also been used for modelling the covariance of lower-order cumulants (Munshi & Valageas 2005). 

The statistics of the projected convergence field can be constructed using a suitable defined variable r) = (k— ftmin)/ Kmin where k(Q, r s ) = 
— f Q e dr w(r,r s ). The variable 77 follows the same statistics as the density parameter S and under some simplifying assumptions and using 
hierarchical ansatz it can be shown that = S^/rj and similar results also hold at higher order i.e. = Kl/rf . The overall dependance on 
the cosmology is absorbed in the definition of rj and the skewness S3 = 3Q, kurtosis K\ = 47i a + 12Rt ~ I6Q4 parameters, defined in terms 
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Gravity induced Non-Gaussianity vs Redshift - Halo Model 
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Figure 5. The Halo Model is used to predict the skew spectra (left panel), Sj (middle panel) and Sj (right panel). The source redshift is unity. The 
underlying background cosmology is that of WMAP7. No smoothing window was assumed. A sharp cutoff at l m ax = 2000 was used for this calculations. Halos in 
the mass range of 10 3 Mq — 10 le A/Q were include in this calculation. The halo model expression for the bispectram is defined in Eq. d29t - Eq.lfJTJ. The dashed 
lines correspond to the analytical model prescribed in Eq.( l341 . 

of the hierarchical amplitudes, Q3 and R a , Rb are insensitive to the background cosmology (Munshi & Jainl200ll ; lMun shi 2000; Munshi & Coles 
l2000h . 

5.2 Halo Model 

The Halo Model relies on a phenomenological model for the clustering of halos and predictions from perturbative calculations on large scales 
to model the non-linear correlation fun ctions. The halo over-de nsity at a given position x, S h (x, M; z) can be related to the underlying density 
contrast <5(x, z) by a Taylor expansion (M o. Jing & Whit e! 1997b . 



<S h (x, M; z) = 61 (M; z)5(x, 2) + ^b 2 (M, 2)5 2 (x, z) + 



(27) 



The expansion coefficients are functions of the threshold v c = S c /a(M, z). Here S c is the threshold for a spherical over-density to collapse and 
cr(M, z) is the rms fluctuation within a top hat filter. The halo model incorporates perturbative aspects of gravitational dynamics by using it to 
model the halo-halo correlation hierarchy; the nonlinear features of this take direct contributions from the halo profile. The total power spectrum 
P l (k) at non-linear scale can be written as 



l!(k,k); P 2h {k) = [ll{k)fP{k); 



P 2h (k)+P lh (k) 



(28) 



(Seliak bOQOl) . The minimum halo mas that we consider in our calculation is 10 3 Mq and the maximum is 10 16 M . more massive halos do 
not contribute significantly owing to their low abundance. The bispectram involves terms from one, two or three halo contributions and the total 
can be written as 



B t {k l M,k 3 ) = B 3h (fc 1 ,fc 2 ,fc 3 ) + B 2 ' 1 (fei,fe2,fc3) + B 1 ' 1 (fci,fc2,fc3); 
B lh = lt{kx,k 2 M)\ 



B 2h (ki,k 2 , k 3 ) = I 2 {ki, k 2 )Ii(k 3 )P(k 3 ) + cyc.perm. 



B 3h (ki,ki,k 3 ) = [2J(k 1 ,k 2 ,k 3 )ll(k 3 ) + I 2 (k 3 )}ll(k 1 )ll(k 2 )P(k 1 )P{k 2 ) + cyc.perm. 



(29) 
(30) 
(31) 



The kernel J(ki, k 2 , k 3 ) is derived using second-order perturbation theory dFrvll984 Bouchet et a j|l992h and he integrals lf L can be expressed in 
terms of the Fourier transform of halo profile (assumed to be an NFW dNavarro, Frenk & Whitell 19961) : 



/„''/■■ 1 ■/•- /•-,.::>= / ( dn< ^ Z) MM)y(ki,M)...v(k^M); y(k, M) = JL j** dr4nr 2 p(r, M) 



sin(fcr) 



kr 



(32) 



© 0000 RAS, MNRAS 000, 000-000 



10 Munshi et al. 



Gravity Induced Non-Gaussianity 
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Figure 6. The skew-spectra and Sl' J are plotted for z s = 1 as a function of wave number I. The underlying cosmology is that of WMAP7. A tophat 

window is assumed. Various curves correspond to different smoothing angular scales as indicated. The resolution is fixed at l m ax = 4000. The smoothing angular 
scales considered are 6 S = 5' (solid lines), 8 S = 25' (long-dashed lines) and 9 a = 55' (short-dashed lines) respectively. The skew-spectra are defined in Eq. )19> . 
The underlying bispectrum is constructed using the analytica model prescribed in Eq. d34t . It is interesting to note that at smaller I the skew-spectra with larger 
smoothing angular scales 8 S dominates. However smaller smoothing angular scales dominates at higher / resulting in a higher values of the corresponding one-point 
skewness parameters. 



;(2) 




Gravity induced Non-Gaussianity vs Redshift 
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cosmology is that of WMAP7. A tophat window is assumed. Various curves correspond to different smoothing angular scales as indicated. The resolution is fixed at 
Imax = 4000. The smoothing angular scale is fixed at 8 S = 25' . Notice that use of broader window not only removed power at smaller angular scale, it also changes 
the overall normalisation of the skew-spectra. The skew-spectra for any specific smoothing angular scales increases with lowering of the source redshift. This is due 
to the fact the PDF of convergence for higher redshift is more Gaussian than at a lower redshift. At a lower redshift the highly evolved large scale structure results in 
higher departure of the convergence statistics from Gaussianity. 



The mass function is assumed to be given by the Press-Schechter form dPress & Schechterll974 . The results are obtained by using Eq.d29ll-Eq.01t 
in Eq.dl2b and Eq.dl3b. The convergence power spectra and bispectra thus computed are then inserted in Eq.dl9j. 

Results from the Halo Model analysis are plotted for the skew-spectra S^ ' (left panel), Sjy (middle panel) and S^ 2 ' (right panel). The 
source redshift is fixed at unity. The underlying background cosmology is that of WMAP7. No smoothing window was assumed. A sharp cutoff at 
Imax = 2000 was used for this calculations. As mentioned, halos in the mass range of 1O 3 M - 1O 16 M were included in this calculation. The 
halo model expression for the bispectrum is defined in Eq . d29 b - Eq.d31l>. 
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Primordial Non-Gaussianity (Local Type) 
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Figure 8. The skew-spectra Sj^ and ai'e plotted for source redshift z s = 1.0 as a function of wave number I. The skew-spectra correspond to the 

primordial bispectrum of local type Eq. d391 . The normlaisation coefficient is set to unity /fe£ = 1. The underlying cosmology is that of WMAP7. A tophat window 
is assumed. Various curves correspond to different smoothing angular scales as indicated. The resolution is fixed at l max = 2000. Three angular scales that we plot 
correspond to 9 S = 5' (solid lines), 9 S = 25' (long-dashed lines) and 9 S = 55' (short-dashed lines) respectively. 
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Figure 9. The skewspectra Si°^ , Sj^ and 5 ; ' 2 ' are plotted. The smoothing angular scale is 9 S = 25'. The three different models that are depicted are local (solid 
lines), equilateral (long dashed lines) and folded (short dashed lines) models. These models are described by Eq. J39h Eq. )4U and Eq.{43) respectively. The redshift 
is fixed at unity z s = 1. All of the non-Gaussianity parameters describing various models are fixed at unity i.e. /«?£ = 1, fjjjf = 1 and fffjf = 1. 



5.3 Perturbative calculations in the quasi-linear regime and their extensions 

In the weakly non-linear regime (S < 1), the description of gravitational clustering can be described by perturbation theory Bernardeau et al I d2002h . 
However, the perturbative treatment breaks down when density contrast at a given length scale becomes nonlinear (S > 1) which significantly 
increases the growth of clustering. Perturbative studies of gravitational clustering have attracted a lot of attention. Starting with lPeeblesI fl980). 
there have been many atte mpts to reproduce the c lustering of a self-gravitating fluid in a cosmological setting is typically tackled by brute force 
using N-body simulations iBernardeau et a il l l2002l) . Expanding the density contrast in a Fourier series, and assuming the density contrast is less 
than unity, for the pertubative series to be convergent, we get 

5(k) = 5( 1) (k) + 5 (2) (k) + 5( 3) (k) + ... ; S( 2 \k) = |^l|^fc(k 1 +k 2 -k)F 2 (fc 1 ,fc 2 )«sW(k 1 )5< 1 »(k 2 ). (33) 

The linearized solution for the density field is (^'(k); higher-order terms yield corrections to this linear solution. Using a fluid approach known 
to be valid at large scales (and before shell crossing) one can write the second order correction to the linearized density field using the kernel 
F2 (ki , k 2 ). Newtonian gravity coupled to the Euler and continuity equation employed to solve a system of non-linear coupled integro-differential 
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Figure 10. The skew-spectra , S^' and S, are plotted for two different source redshifts z s = 1.5 (long-dashed lines) and z s =0.5 (solid-lines) as a function 
of the wave number I. The lower pair of curves in each panel correspond to primordial bispectrum of local type. The gravity induced skew-spectra that dominate all 
scales at each redshifts are also shown in comparison (upper pairs of curves). 



equation reproduces the kernels F2(fci, fc 2 ) Fs(ki, &2, ^3) when solved perturbatively order by order. The expression for the matter bispectrum 
can be written in terms of an effective fitting formula that can interpolate between quasilinear regime and the highly nonlinear regime: 



Bi(ki,k 2 ,k 3 ) = 2F 2 (ki,k 2 )Pz'L(ki)Pi'L(k 2 ) + cyc.perm.; 

F 2 (ki, k 2 ) = ja(n e , k)a(n e , k) + ^~~7^|~~ + ^J^ J b ( Ue > k ) b ( n c k ) + f y lifa 



c(n e , fc)c(n e , k) 



(34) 



The coefficients a(n e , fc), 6(n e , k) and c(n e , k) are defined as follows: 



a(n e , k) = 



1+^8 



1 + (g/4)» 



:+3.5 



b(n e , k) = 



l + 0.4(n e + 3)g" e+3 

1 _|_ q n ef f+3 



;(n e , k) = 



1 + (2?) 



Be+3.1 



1 + 



4.5 



1.5 + (n e + 3) 4 



(35) 



jd In k, q is the ratio of a given length scale to 
2"")/{l + 2 n "). Similarly a 8 (z) = D(z)a 8 . 



Here n e is the effective spectral slope associated with the linear power spectra n e = d In P l 
the non-linear length scale q — k/k n i, where k 3 /2iv 2 D 2 (z) Py ln (k n i) — 1 and Q 3 (n e ) = (4 
At length scales where 1 which means the relevant length scales are well within the quasilinear regime a — b — c = 1 and we recover the 
tree-level perturbative results. In the regime when g > 1 and the length scales we are considering are well within the nonlinear scale we recover 
with b 



2 (z)^0.7Q 3 (n e 



c — 0. In this limit the bispectrum becomes independent of configuration and we recover the hierarchical 



form of bispectrum discussed before. However is there are weak violations of hierarchical ansatz in the highly nonlinear regime is still not clear 
and can only be determined higher resolution N-body simulations when they are available. Similar fitting functions for dark energy dominated 
univer se calibrated aga inst simulations are also available and at least in the quasilinear regime most of the difference comes from the linear growth 
factor dMa et al]|l999l) . The analytical modeling of the matter bispectrum presented here is equivalent to the so-called Halo Model predictions 
presented above. 

We have used this model to construct the analytical predictions for various skewness parameters and the corresponding skew-spectra. The 
results are plotted in Figures [3] [4] [6] and [7] In Figure [3] we have plotted the three skewness parameters 5 1 ' ', S 1 ' 1 ' and as a function of the 
smoothing angular scales 8 S as defined in Eq.(|9}. In Figure-|4]we change the source redshift to compare predictions. In total we compare three 
different redshifts z s — 0.5, 1.0 and z a — 1.5 respectively. We use top hat filters with different angular smoothing scales. The skew-spectra, 
defined in Eq.(|19l>, are integrated measures and their value at a specific harmonic depends on modelling of the bispectrum for the entire range of 
harmonics being considered. The skew spectra are plotted as functions of harmonic / in Figure|6]and Figure|7]respectively. In Figure|6]we consider 
the redshift z s — 1.0 and in Figure[7]results for two different redshifts, z 3 = 0.5 and z 3 = 1.5, are compared for a given angular smoothing. The 
oscillatory behaviour seen in these plots is due to our choice of liter function i.e. top hat window. 
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5.4 Primordial non-Gaussianity: bispectrum 

A recent (controversial) claim of a detection of non-Gaussianity dYadav & Wand elt 2008) in the Wilkinson Microwave Anisotropy Probe 5-year 
(WMAP5) sky maps, has boosted interest in cosmological non-Gaussianity. Much of the interest in primordial non-Gaussianity has focussed 
on a phenomenolo gical 'local /jvx' parametrization in terms of the perturbative non-linear coupling in the primordial curvature perturbation 
dVerde et alj2007b : 

= $l{x) + (a;) - <*£(*)», (36) 

where &l,(x) denotes the linear Gaussian part of the Bardeen curvature and /jvl is the non-linear coupling parameter. A number of models have 
non-Gaussianity which can be approximated by this form. The leading order non-Gaussianity present in this model is at the level of the bispectrum, 
or in configuration space at the three-point level. M any studies in volving primordial non-Gaussianity have us ed the bispectrum, motivated by the 
fact that it contains all the information about Jnl jB abich 2005!). This model has been extensively studied dKom atsu, Sp ergel & Wandeltll2005l : 
ICreminellil2003l ; ICreminelli et al.l2006l;,Medeiros & Contaldil2006l : ICabella et al. l2006l ; ISmifh. Senatore & Zaldarriagal2009l) . with most of these 
measurements providing convolved estimates of the bispectrum. It is interesting to note here in the context of bispectrum estimation from CMB 
sky, optimized 3 -point estimators were introduced by Heavens ( 1998), and have been successively developed dKomatsu. Spergel~& Wandelt [2005l : 
ICreminelli et al J 120061 : ICreminelli. Senatore. & Zaldarriaga 2007; S mith. Zahn & Por e 2000; Smith & Zaldarriaga 2006) to the point where an 
estim ator for Jnl which saturates the Cramer-Rao bound exists for partial sky coverage and inhomogeneous noise dSmifh. Senatore & Z aldarriaga 
2009). Approxi mate forms also exist for equilateral non-Gaussianity, whi ch may arise in models with non-minimal Lagrangian with higher- 
derivative terms dChen. Huang & Kachrul 12006. ; Ich en. Easther & Li ml l2007l) . In these models, the largest signal comes from spherical harmonic 
modes with £\ ~ £2 — £3, whereas for the local model, the signal is highest when one £ is much smaller than the other two - the so-called squeezed 
configuration. 

In the Fourier space the primordial bispectrum of local type defined in Eq.d36l> takes the following form: 



Bi oc (kiki,k3) = VnI [i&(fci)i&(fca) + cyc.perm.] . (37) 

The primordial potential power spectrum in standard inflationary models takes a power law form P*(fc) oc k"~ 4 . In t he linear regime, the 
primordial bispectrum for the density field B]? rlm in case of local model evolves according to the following expression (see lHikage et al.l (2006) 
for a detailed derivation and discussion): 



Bj oc (kiki, k 3 ; z) 



r) r\oc 

D(z) 



-Piin{ki,z)P lin {k 2 ,z) + cyc.perm. 



M{k 1 )M(k 2 ) 



S(k,z) = D(z)M(k)$Qt,z); M(k) = - 



k T(k) 

rn. 



(38) 
(39) 



Here D(z) is the li near growth factor normalised such that D(z) — > 1/(1 + z) and T(k) is the transfer function given by an approximate 
expression found in BBKsj dl986, ). According to standard inflationary predictions P* (k) oc k n ~ A and the linear power spectra for the density 
is given by Pf in {k, z) = D 2 (z)M(k) 2 (k) . The primordial bispectrum for the density can similarly be expressed in terms of that of the 
primordial potential perturbations Bg" m (k i : k 2 , fc3, z) = D {z) M(k\)M(k2)M(k-j)B 1 ^ zm (ki, k 2 , ks). The primordial potential bispectrum 
for the equilateral type can be expressed as (Creminel li et al J 20061) : 



Bl qul = 6f N T [-(P*(fa)P*(fe) + cyc.perm.) - 2(P $ (fc 1 )P $ (fc 2 )P*(fc 3 )) 2/3 + (Pj /3 (fci)P| /3 (fc 2 )P*(fc 2 ) + cyc.perm.)} 



(40) 



The primordial density bispectrum for the equilateral case B^ qu% can be expressed, following the same procedure that we followed for the 
local type, as: 



-2 (M(fci)A4(fe 2 )M(fc 3 ))^ 1/3 {PL(ki,z)PL(k 2 ,z)PL(k 3 ,z)} 2/3 

+ \ ^StTT, I {PL{k 1 ,z)PL{k 2 , zfPhfa, Z ) 3 } 1/3 + cyc.perm) } . (41) 
\M(k 2 ) 1/3 M(k 3 ) J J 

However in contrast to the local model , in the equilateral model the functional form of the expressions do not have any connection to fundamental 
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physics but are just fits where the exact expressions are more complicated. The folded or flattened model that is maximized for k% w £3 w fei /2 
is well approximated by the following form: 



B J 



" 6/j(/f [(P*(fc 1 )P*(fc 2 ) + cyc.perm.) + 3{P^(k 1 )P^(k 2 )P^(k 3 )) 2/3 - (P^ /3 (fc 1 )P| /3 (fe 2 )P 4 ,(fc 2 ) + cyc.perm.)] (42) 
In terms of the density perturbations, we get the following expression: 



■ s fk!k! - k "->- D(z)i{M(k 1 )M(k 2 ) 



[ ( X (ki)M(k 2 ) P '' n (fcl ' z) P "" (fe ' z) + c ^ erm -) 

3(M(fc 1 )^(fc2)^(fc 3 ))- 1/3 {PL(fel^)PL(fo,z)PL(fc3^)} 2/3 



A<(fc0 1/3 {p f w(fc z)P f in {k 2t z) 2 PL(k 3 , zf} 1/3 + cyc.perm. ) | . (-13) 

A1(fc 2 ) 1/;! X(fc3) 

The folded or flattened form of bispectrum appears in canonical single field models where the initial Bunch-Davies vacuum is modified. 

The evolution of the primordial bispectrum is different to that generated by gravitational evolution. The angular dependence for the 
gravitationally-induced bispectrum is also different. On large angular scales, which will be probed by future weak lensing surveys, gravitational 
instability may not have erased the memory of primordial non-Gaussianity, which can provide supplementary information to results obtained from 
CMB observations. 



6 ESTIMATORS AND THEIR SCATTER 

As noted above, the estimators for the skew-spectra can be most easily computed by cross-correlating maps in the harmonic domain. These maps 
are constructed in real space by applying various derivative operators. The recovered skew-spectra will depend on the mask, if one is present, 
because a mask typical ly introduces mode-mode coupling. The approach we adopt here to reconstruct the unbiased power spectra in such a case is 
the Pseudo-C; method jHivon et alj|2002h . This approach depends on expressing the observed power spectra Ci in the presence of mask as a linear 
combination of unbiased all-sky power spectra. 

The three different generalized skew spectra that we have introduced here can be thought as cross-spectra of relevant fields. We denote these 
generic fields by A and B and will denote the generic skew-spectra as s\ A,B \ The skew-spectra recovered in the presence of masks will be 
represented as s\ A ' B ^ and the unbiased estimator will be denoted s\ A ' B K The skew-spectra recovered in the presence of mask s\ A,B ^ will be 
biased. However to construct an unbiased estimator s\ A , B ^ for the skew-spectra the following procedure is sufficient. The derivation follows the 
same arguments as detailed in Munshi .~Smidt & Cooravl j2010h and will not be reproduced here. 

S lA ' B] = -±-J2 Al ™ZL; §\ A ' B] =Y,M W S\ A ' B] ; M w = ^—J^Iun-M 2 ; {A,B}e{K,K 2 ,(X7K-X7K),X7 2 K}. (44) 

m I' + I'l" 

The mode-mode coupling matrix M is constructed from the power spectra of the mask w^i and used for estimation of unbiased skew-spectra S A ' B . 
Typically the mask consists of bright stars and saturated spikes where no lensing measurements can be performed. The results that we present here 
are generic. The estimator thus constructed is an unbiased estimator. The computation of the scatter covariance of the estimates can be computed 
using analytical methods, thereby avoiding the need of expensive Monte-Carlo simulations. The scatter or covariance of the unbiased estimates 
(5S A,B 8S A ' B ) is related to that of the direct estimates {5S A ' B 5S A ' B ) from the masked sky by a similarity transformation. The transformation is 
given by the same mode coupling matrix M: 

SI A ' B] = ^[M-V^' 51 ; {SS lA ' B] SS\f' B] ) =Y / Mr L \6S [ ^ B] 6S [A ,- B] )Mj l ,-, {S lA ' B] ) = S lA ' B] . 

I' LL' 

The power-spectra associated with the MFs are linear combinations of the skew-spectra (see Eq.([3}). In our approach the power spectra associated 
with the MFs are secondary and can be constructed using the skew-spectra that are estimated directly from the data. 

No construction of an estimator is complete without an estimate of its variance. The variance or the scatter in certain situations can be 
computed using Monte-Carlo (MC) simulations which are computationally expensive. In our approach, it is possible to compute the covariance of 
our estimates of various Sis, i.e. {SSiSSf) under the same simplifying assumptions that higher-order correlation functions can be approximated 
as Gaussian. This allows us to express the error covariance in terms of the relevant power spectra. The generic expression can be written 
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[Vt 2 \ = Et M ~V I^ (2) ]«; <^f W') = E WL{6[Vr]i6{V^)i>}M-,\, (45) 

We would like to point out here that, in case of limited sky coverage, it may not be possible to estimate the skew-spectra mode by mode as the 
mode coupling matrix may become singular and a broad binning of the specra may be required. 

({8S\ x ' Y] ]SS l x- Y] }} = f^^i [c\ X ' X] C l v' Y] + [SI X - Y] ] 2 } 6 U >; {X, Y} e {k, k\ Vk(Q) ■ Vk(A), (46) 

(k 2 k) (k 2 V 2 /i) 

Here the fraction of sky covered by the survey is denoted by / s ky. The expressions for the skew-spectra are quoted in ' and 

D; are given in Eq.GUl. The expressions for covariance also depend on a set of power spectra i.e. S{ ' , £/ ' , 

and S*' K . These are given by the following expression: 

c v-v,v.v = y^c v + AT,,)(C,„ + JV,«)[ii(ii + 1) + la(Ia + 1) - l(' + l)?Ju>i»W v W vl \ (47) 
i'i" 

C, [k2 ' k2] = E^' + + N v ,)ll, v ,W v W v r, C |v 2 -.v 2 K ] = Z 2 (z + 1)2(C; + (4g) 

t'i" 

Here Ci is the ordinary "theoretical" convergence power spectrum defined in Eq. d 1 2t that includes noise i.e. Ci Wi is replaced with C{W\ + Ni with 
TV; = 4.7rof /Ng a i. Here ai is the intrinsic ellipticity distribution of galaxies and N ga i is the number of galaxies per arc-minute square. Using these 
equations it is possible to compute the scatter in various skew-spectra. These results can also be extended to take into account the cross-correlation 
among various skew-spectra extracted from the same data. Using a less compact notation we can write 

(Ssf^Ssf^) = / s - ky ^ [cf^C^ + [sf^f] (49) 

(ssl^ssf^) = f£^L- [ C [«v 2 ] c rv.v,v.v, + [5 |^,v^ ]]2 ] (50) 

{(55 [v K .v K ,v 2 K ] 5s [v K .v K ,v 2 K]) = J c [v 2 K ,v 2 .] c [v.v,v.v] + [s [v K .v.,v 2 K]]2 j _ (51) 

The cumulative signal to noise upto a given l max using these expression for estimators can now be expressed as: 



N 



= f sky g (21 + 1) { - ( f' } 2 1 (52) 

The signal-to-noise for the other two estimators S^ 1 ' and 5 1 ' 2 ' can be defined likewise. The different skew-spectra that we have studied here are 
not completely independent. Their covariance can be analysed using the same procedure, allowing their joint estimation from a single data set. 

(Ssf^Ssf^) = /- ^ [ C f^c\^ + sf^s\^] (53) 

<^ ^ lK2 ' Kl ^S ^ [VK ' VK ' v2 ' t, } = ./Tw^— f S [«^ 2 ^M».v B i + ( ^,v K .v.) (B ,v^)l (54) 

y 21 + 1 L J 



(55) 



The above results can be generalized to compute the cross-covariance of Si from different sources of bispectrum. The following quantities are 
required to compute the necessary cross-covariances. 

c [«,v 2 *] = + 1)Ci . c [ K 2 ,v K .v K] = ^(CvWv + Ni){Ci„Wv + Ni„)lln»[l\l' + 1) + + 1) - 1(1 + 1)]; (56) 
C [»,v«.v»] = _ ^p# (J , + 1} + + l{l + i )]Bwi „j m „ Wl , Wv , (57) 

VI" 
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Figure 11. The cumulative signal S/N(> I) to noise associated with gravity induced skew-spectra for Sj°\ and S?P^ ai'e plotted as a function of the wave 
number I. We have assumed a full sky coverage f s ky = 1. The results plotted are for l m ax = 4000 and the smoothing angular scale is 9 S = 5'. The curves from 
top to bottom are Sj , and respectively. The signal-to-noise is highest for fiW. The non-Gaussianity decreases with the increase in source redshift. 
However the power spectrum and hence the scatter increases with redshift. This makes is easier to probe non-Gaussianity at relatively lower redshifts. 



We have discussed the lowest-order departure from Gaussianity in MFs using a third order statistic, namely the bispectrum. The next-to- 
leading descriptions are characterized by the trispectrum which is a fourth order statistics. It is possible to estend the definition of skew-spectra 
to the case of kurt-spectra or the power spectrum associated with tri-spectra. The power spectra associated with the Minkowski Functionals can 
be defined completely up to fourth order using the skew- and the kurt-spectra. However, the corrections to leading order statistics from kurt- 
spectra are sub-dominant and leading order terms are consequently sufficient to study the departure from Gaussianity. In any case it is nevertheless 
straightforward to implement an estimator which which will estimate the power-spec trum associated with the MFs f rom noisy data by including 
both third order and fourth order statistics; this issue has been dealt with in detail in (Munshi, S midt & CooravllioiOl) in the context of CMB sky. 
The same results will also be applicable for weak lensing surveys. 

In addition to the three genralised skew-spectra that define the MFs at lowest order in non-Gaussianity, it is indeed possible to construct 
additional skew-spectra that work with different set of weights. In principle arbitrary number of such skew-specta can be constrcuted though they 
will not have direct links with the morphological properties that we have focussed on, in this paper they can still be used as a source of independent 
information on the bispectrum and can be used in principle to separate sources of non-Gaussianity, whether primordial or gravity induced. 



7 NEXT TO LEADING ORDER CORRECTIONS TO THE MINKOWSKI FUNCTIONALS FROM THE TRISPECTRUM 

The skew-spectra Sj and the related skewness parameters completely specify the MFs at leading order. The next-to-leading order corrections 
are determined by a set of four kurtosis parameters K^' . These generalised kurtosis parameters are constructed from the trispectrum using varying 
weights to sample different modes. This method is very similar to construction of the generalised skew-spectra and their associated skewness 
parameters from the bispectrum described in previous sections. The kurt-spectra are constructed by cross-correlating maps that are constructed 
from original maps and combinations of maps constructed from the original map e.g. V$(S1) and V 2< 3?(f2). The four kurtosis parameters are 
natural generalizations of the ordinary kurtosis and can be most easily be esti mated in real spac e. The normalization of these kurtosis 

parameters are determined by suitable combinations of powers of parameters <jq and ai dMatsubaraluOlOf) . 



K (o )= i _ (* 4 (n))c , K m_ i g (»»v»») = (4fj)vx% , (58) 
K (2) _ K (2a) + K m s i g («iv»pv»») + _J_^(|vk|4) = M(v«("))l 2 (vM")))c + (l(v«)l 4 ) c . (59) 

Jf (3) = ;-j-^if (|VK|4) = <|V f ( "y )c ; where \X7k(£i)\ 2 = Vfs(A) ■ V«(fi). (60) 

ACT q (T-^ ACT qCT i 

Unlike the skewness parameters the kurtosis parameters get contributions also from Gaussian (unconnected) components. The subscript c above 
however refers to the non-Gaussian or the connected part of the contribution which is directly liked to the trispectrum. 
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The correction to the Minkowski Functionals SV^ (y) as defined in Eq.lO from the n ext to leading order terms consists of both the Kurtosis 
parameters as well as the product of two skewness parameters S (i) jMatsubarj|20ld) : 



72 

[g (0) ] 2 
72 

[gW] 2 
72 



H (y) + 
H r (v) + 



24 

K (o) _ s (o) s (i) 



24 

K W _ 5 (0) 5 (i) 



24 



3 r ^(i) l2 



H 2 (v) 



K 



(i) 



(0) o(2) 



K 



(2) 



(l)c(2) 



(61) 



(62) 



The analytical modelling of four-point correlation functions is most naturally done in the harmonic domain. They are described by the angular 
trispectrum T^{L) which is defined through the relation (K iimi K i2m2 Ki 3m3 K !4m4 ) c = J2l hihL Ii 3 UlT '^ (L ). The trispectru m 7^ 2 (I) 
is expressed in terms of the reduced trispectrum f/ 1 / 2 (£). Following expression was introduced by dHull2000l 2001| : IHu & Okamotoll2002f ) and 
encodes all possible inherent symmetries. 



The matrices in curly brackets represent 6j symbols which are defined using 3j symbols; see lEdmondsl Jl968h for more detailed discussions. 
The entities P^i^ (0 can be further decomposed in terms of the reduced trispectrum (I). A specific model for the non-Gaussianity - either 
primordial or gravity-induced - has a specific prescription for the reduced trispectrum which in turn describes the next -to-leading-order corrections 
to the MFs. 



\h+h 



h 
h 



h 



h 

la 



h 
U 



(63) 



Pllll (0 = AllHl) + (-lf u T%t(l) + (-l) Et T^(0 + (-lf^v T hh (l) . x L=h+h + ,. X u= l, + l 4 + l. (64) 

In addition to the original convergence trispectra [T' '] generally used in the literature we can define a set of four trispectra which uses different 
weights to samples of modes deefined by the quadruplet of harmonic numbers U . 

[r (o) ]^:(/) = it 3 \ i :(iy, [T (i) ]^:m = \ i^i + 1) + uk + i) + hi 3 + 1) + hu + iw^ny, (65) 

[T {2) ]\\\\{1) = \ [1(1 + 1) - (h(h + 1) + k(k + l))(l a (/a + 1) + h{h + 1))] 7^(0; (66) 

[7" (3) ]^:(0 = \ {(h(h + 1) + l a (i a + 1) - 1(1 + + 1) + I*(f4 + 1) - 20 + 1))] 71^(0- (67) 

The power spectrum associated with these kurtosis parameters, i.e. the kurt-spectrum, is defined in terms of the trispectrum, extending the previ- 
ously defined skew-spectra along fairly obvious lines. The estimation of these kurt-spectra would be done by cross-correlating relevant fields used 
to construct the related kurtosis in real space : 

^'^EE^^^wWi =Y}^]\±{l)Ji ll2l Ji 3 i^ ^(2Z + l)if«=Jf«. (68) 

The error and covariance associated with these kurt-spectra can be computed using exactly the same formalism that we described in the context of 
estimation of skew-spectra. The kurt-spectra being a power spectra will contain more information compared to the kurtosis which is a one-point 
estimator. Though one may be interested in principle to extract the entire trispectrum, it may be more realistic to use the kurt-spectra because of 
the likely low signal-to-noise associated with individual harmonic modes. 

To compute the kurtosis one needs a resonable model to compute the trispectra for the convergence field. This is typically done using the 
paraphernalia of the Halo Model we introduced before. The modelling of the gravity-induced trispectrum in the Halo Model follows the same 
principle as before. It involves contributions from one-, two-, three- and four-halo contributions and the total can be written as: 

T s (k 1 ,k 2 ,k a ,k 4 ) = Tl h (k u k 2 ,k 3 ,k 4 ) + Ti h {k 1 ,k 2 ,k 3 ,k 4 )+Ti h (k 1 ,k 2 ,k a ,^ (69) 
The expressions for various contributions are listed below. These can be expressed in terms of k 2 , . . . , fc M ; z) defined above. Notice that 



© 0000 RAS, MNRAS 000, 000-000 



1 8 Munshi et al. 



z=0.5 



z=1.0 





6 


A 


4 








2 








i 


1111 






/ 

— / 


6 = 5 Z 

s 


/ 




/ 




_/ 




-I 
t 




I 1 ! ! I ! ! 1 1 



1000 2000 
I 




1000 2000 
I 



15 
10 



z=1.5 



> 5 



_ / 
/ 

/ 

_ / 



~~\ i j i i 

1000 2000 
I 



Figure 12. Same as the previous figure but galaxy shot noise is included in the computation of scatter. 

two-halo (2h) term has two contributions. In the case of one represented by T 3 i there are three points in the first halo and one in the second. 
Whereas T 22 represents two points in each halo. 



Tg = I°(ki,ki, k 3 , fc 4 ); 



Tii = Pu„(ki)ll(k2, k 3 , fc 4 )/i 1 (fci) + eye. 



perm. 



(70) 

T 2 2 2 = PL{ki2)ll(ki,k2)ll(k 3 , k 4 ) + cyc.perm.(71) 



Tg h = Bi in (k 1 ,k 2 ,k 3 )l2(k 3 ,k4)ll(k 1 )ll(k 2 ) + Pun(ki)Pu„(k2)l'i(k 3 ,k i )I 1 (k 1 )ll(k 2 ) + cyc.perm. 
Tt = ll(k 1 )ll(k2)ll(k 3 )ll(k 4 )TL(ki,k 2 ,k 3 , kt) + ll(k 4 )ll(k 1 )ll(k 2 )ll(k 3 ) + cyc.perm 



(72) 
(73) 

Here Pf in is the linear power spectrum for the density contrast 8 and Bf in (ki, k 2 , k 3 ) and Tf in (ki, k2, k 3 , k 4 ) is the tree level expression for 
the bispectrum in quasilinear regime. The general express ion for J M (fci, k2, ■ ■ ■ , fc^; z) is quoted in the expression Eq.((32}- Detailed derivations 
and discussions of these expressions can be found in e.g. ICoora vl fcOOlk lc oorav & Sethi II2002I) . Extension of perturbative approaches can also 
be employed for computation of gravity induced trispectrum. The accuracy of any analytical modeling is more difficult for the higher order 
multispectra and depends by and large on more inputs from numerical simulations. The projected tripsectrum or the convergence trispectrum 
nail (0 can b e ex P resse d in terms of the underlying mass trispectrum (I): 



hall (!) = IhhllhUl 



d\{r) 



h 



I2 



u 



d,A (r ) ' d,A (r) ' gU (r ) ' <1a (r ) 



(74) 



The trispectrum for primordial non-Gaussianity exists in the literature for the local model . The results presented here clearly are generic and can 
be deployed to analyze arbitrary models. It is worth mentioning here that while modelling of trispectrum is relevant for computation of corrections 
to the leading order terms they are also important in modelling the scatter in computation of ordinary power spectrum. Hence the errors in ao and 
(Ji e.g. will involve the one-point kurtosis parameters K'^ if contributions from non-Gaussianity are taken into account. The correlation functions 
that represent these kurt-spectra in real space are constructed using derivative operators on the original convergence map and can be useful for 
surveys with smaller sky coverage. 



# w (ni,n 2 ) = (k'(cii))k(£i2))c; K m (Qi,n 2 ) = (^(fiOMf^v^na)]) 



A: (2) (ai,o 2 ) = (v«;(fii) ■ Vk^ok^v^^ODc x w (ni,f2 2 ) = (V/t(fii) • VK(fii)[VK(fi 2 ) • v«(n 2 )]> 



'(3)/ 



(75) 
(76) 



These correlation functions can be computed directly in real space without any harmonic decomposition. 



8 CONCLUSION 

Weak lensing observations offer the potential to probe the cosmological density distribution in an unbiased way. Since the angular scales probed 
by weak lensing are sensitive to non-Gaussianity, primarily that generated by gravitational clustering, this technique offers us the chance to push 
our understanding of the statistical properties of the cosmological matter field far beyond current limits. 
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The statistical characterization of gravitational clustering is most often performed using a hierarchy of higher order correlation functions or 
their collapsed counterparts which correspond to the moments of the convergence field k. However, it is well known that non-Gaussianity can 
also modify the morphological properties characterized by the MFs of the relevant field k. The MFs therefore encode information about the non- 
Gaussianity and can be used as an estimator. At leading order the MFs depend on three generalized skewness parameters S°, S 1 and S 2 . These 
parameters are one-point statistics constructed from the bispectrum -B; 1 ; 2 ! 3 using different weights for individual modes. We have generalized 
these one-point estimators to a set of power spectra, namely Sf , and Sf . We studied how they can be expressed in terms of the bispectrum 
Bi^^. In real space these power spectra are related to the relevant correlation function Eq. J23t . Though the correlation functions associated with 
the skew-spectra are two-point statistics in terms, of spatial order, they actually are third (lowest) order in terms of non-Gaussianity. Hence they 
carry information about the bispectrum. These statistics are in fact known as the cumulant correlators and the first of these statistics, S^°\ is 
already well studied in the literature. The expression for a generic cumulant correlators of order p + q is (k p (S7i W^fj^Vlt probes multispectra 
of order p + q and are known to be related with bias associated with over dense objects in 3D or hot-spots in 2D ( Munshi 2000). 

The skewness parameters define the leading-order terms to the MFs. The next-to-leading-order terms are associated with the convergence 
trispectrum. The convergence trispectrum in turn is expressed in terms of trispectrum of the projected density field. The generalized kurtosis 
parameters and their related power spectra can likewise be constructed from the convergence trispectra. The corresponding representations in the 
Fourier domain are named as the kurt-spectra. We have not considered these kurt-spectra in our analysis as they are sub-dominant but they can be 
taken into account using the same formalism if required. 

We have shown that the MFs can be decomposed into three different power spectra and that these power spectra can be constructed from 
an equal number of skew-spectra that carry information completely equivalent to the original MFs at the lowest order. These power spectra in 
real space will correspond to correlation functions of fields that are constructed from products of various derivative fields. These spatial derivative 
fields are in turn constructed from the original convergence maps k(8 3 ). These generalized skew-spectra are therefore related to the generalized 
cumulant correlators defined in real space. Each of these skewness parameters can be constructed from the relevant skew- spectra. However, the 
skew-spectra have the greater power in distinguishing different sources of non-Gaussianity. This is related to the fact that individual sources of 
non-Gaussiantity will lead to specific shapes for the skew-spectra that can be tested against the observed data. We have shown that recovery of 
these skew-spectra is relatively straightforward from noisy data and in the presence of a mask. The scatter in these statistics can be estimated under 
certain simplifying approximations. 

In this paper we have initiated a systematic study of these skew-spectra in the context of weak lensing surveys. We have studied how the 
skew-spectra depend on specific choices of non-linearity that include gravity induced non-Gaussianity or primordial non-Gaussianity. We have also 
pointed out that the departure of MFs from Gaussianity is determined by the generalised skew-spectra which are largley independent of cosmology 
but which depend primarily on specific models of primordial non-Gaussianity. The overall amplitudes are determined by the background cosmology 
as they are determined by the power spectrum of convergence. Such a clear distinction promises to help enormously separating the non-gaussianity 
independent of cosmology. 

The formalism we have developed here for the study of non-Gaussianity depends on the well known pseudo Ci approach for the power 
spectrum estimation. In this approach, the effect of any mask and noise can be dealt with in a natural manner. This is achieved using a matrix that 
encodes the mode-mode coupling. We generalized this approach to the context of generalized skew-spectra and showed that the error and their 
covariance can also be constructed in this approach. We also performed a detailed analysis of error characteristics. The analytical characterization 
of errors means numerical costly Monte-Carlo simulations are no longer needed and is a further strength of this approach. 

It is also worth pointing out that although we have considered three generalized skew-spectra which are related to the MFs, it is clearly the 
case that infinitely many such generalized skew-spectra can constructed with arbitrary associated weights that are not directly related to MFs. 
However these generalized skew-spectra can be analyzed jointly to maximize the extraction of the information content. 

One fly in the ointment is that we do not have a complete analytical picture of gravitational clustering. However, a number of variants of 
perturbative techniques which also rely on inputs from numerical simulations are widely in use. We are also reasonably confident that the Halo 
Model is capable of capturing basic features of gravitational instability. We have used these approximations to construct correspondong theoretical 
predictions for the skew-spectra. We study them as a function of redshift of sources as well as the smoothing function to check how sensitive the 
results are to various assumptions about the input physics. 

Non-Gaussianity induced by gravity may be the primary source of non-Gaussianity for weak lensing probes, but recent CMB studies have also 
pointed to the possibility of non-zero primordial non Gaussianity. It is well accepted that CMB studies may be the cleanest probes to primordial 
non-Gaussianity. Nevertheless, large-scale structure probes are known to reach comparable accuracy. It is therefore interesting to see if weak 
lensing observations too can be used to detect and study various models of primordial non-Gaussianity. Motivated by the idea that the skew-spectra 
might be valuable in this direction, we have studied to what extent the skew-spectra can provide valuable information about various models of 
primordial non-Gaussianity. We specifically studied two different models of primordial non-Gaussianity, namely the local model and the equilateral 
models of non-Gaussianity, and compared their contributions against the gravity-induced non-Gaussianity generated due to subsequent evolution 
as a function of redshift as well as angular harmonics. 

The window function which we have considered here is top hat window. Clearly the results can be generalized to any other windows e.g. the 
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Map or Gaussian window functions that too are often used in various observational situations. However the use of different window function is 
not expected to change the overall conclusions. 

We have ignored noise in weak lensing surveys that arises from the intrinsic distribution of galaxy ellipticities. It is expected that noise arising 
from this will somewhat dilute the signatures from the non-Gaussianity, because if increased scatter. However, for a reasonable number-density of 
galaxies, the noise power spectrum will overtake the convergence power spectrum beyond a harmonic mode I where saturation in signal-to-noise 
has already been reached and so will not likely to change the saturation value of the cumulative signal-to-noise. This is true for all of the estimators 
probed as they reach saturation for roughly the same value of I as shown in FigJT] 

Weak lensing statistics are very sensitive to the cut-off in halo mass used in the calculations. We have used halos in the mass range of 1O 3 M 
- 10 16 Mq. Higher-order statistics are typically determined by the high-end tail of the density distribution, i.e. by regions within high mass halos. 
Selective choice of a specific mass range will clearly change the detailed result and can be incorporated in our analysis. The three different skew- 
spectra that we have proposed can be used to separate up to three different components of the non-Gaussianity. Additional skew-spectra can be 
constructed which can be used for a consistency check though they may not have any direct link to the MFs. The results presented here are also 
for a single source plane, e.g. z a — 1, but a realistic redshift distribution of sources can easily be incorporated in our analysis. 

To summarize, we find that for /ml — 1 which specifies the primordial non-Gaussianity the skew-spectrum is typically two orders of 
magnitude lower than the gravity induced non-Gaussianity. This is true for all three different models of primordial non-Gaussianity that we have 
probed irrespective of the source redshift. This will mean for a reasonable value of /ml (say e.g. Jnl~ 100) the gravity induced non-Gaussianity 
and primordial non-Gaussianity will make nearly equal contributions to the various skew-spectra with roughly equal signal-to-noise. The scatter 
does not depend on the model of non-Gaussianity and depend only on the power spectrum. Of the three skew-spectra studied we found that the 
highest signal-to-noise is achieved by the skew-spectra followed by S^°K For all three redshifts we have probed we found that has the 
lowest signal-to-noise and may not be detectable even with all-sky coverage. 

It is worth mentioning here finally that although we have studied the projected or 2D morphology of large-scale structure as probed by weak 
lensing surveys, it is indeed possible to extend these results to 3D weak lensing surveys. The 3D weak lensing survey generalizes the tomographical 
studies to 3D using photometric redshifts. In future many 3D weak lensing surveys will provide us with an unbiased picture of the dark matter 
distribution. Statistical descriptors will be important to quantify such 3D distribution of dark matter. The 3D morpholog y of the large scale structure 
has also been studied extensively studied using morphological descriptors applied to redshift surveys; see (Seth 2006) and the references therein. 
The 3D morphology is far richer than the 2D descriptors considered here for the projected surveys. In 3D there are four MFs which correspond to 
the surface area Vo, volume Vi, extrinsic V3 and intrinsic curvatures V4 respectively. These MFs are used to define various statistics that are linked 
to genus and percolation statistics. Shape statistics have also been introduced to link the MFs with statistical analysis of shapes and are now widely 
used for analyzing galaxy surveys and N-body simulations. 

In future, use of photometric information of galaxies will allow mapping out the dark matter distribution using weak lensing surveys. Such 
3D weak lensing surveys will provide us 3D maps of the dark matter distribution that can be probed using morphological descriptors. The direct 
link with bi- and tri-spectra based approach developed here can be useful in studying growth of structure under gravitational instability. These has 
the potential to greatly enhance the information gained by studying projected catalogs that we have presented here. 
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